clear all
clc

load('D:\My Documents\CIRANO\db_mySQL\tp_2010Dec11_MR7.mat')
conn = database('BDC','','');

yr = 11;
T = 1; % time horizon
R = 0.45; % recovery rate
numSimulate = 300;
%loss_temp = zeros(length(exposureInd)*numSimulate,length(conditions));
loss_portflio = zeros(numSimulate,length(conditions));

for grp = 2: 2 %length(conditions)
    
    condition = char(conditions(grp));
    % calcuate exposure of assigned year
    sqlquery1 = strcat(['select grpRate(modiRate) as grpRate, avg(PPDOS_sum), count(customer) from helena20111207.final_v4_tbl as t where Year(modiFisDate) =  ''' num2str(year_list(yr)) ''' and Month(modiFisDate) = 12 and ' condition ' group by grpRate;']);
    sqlquery2 = strcat(['select customer, grpRate(modiRate) as grpRate, PPDOS_sum from helena20111207.final_v4_tbl where Year(modiFisDate) =  ''' num2str(year_list(yr)) ''' and Month(modiFisDate) = 12 and ' condition ' and modiRate <> -1 group by customer;']);
    exposureAvg = fetch(conn, sqlquery1);
    exposureInd = fetch(conn, sqlquery2);
    
    % set portfolio distribution of rating and exposure
    portfolio (:,:) = cell2mat(exposureAvg(2:numEdges, :));
    display(['Portfolio total firm# in group ' num2str(grp) ' is ' num2str(length(exposureInd))])
    
    
    % get generator for grp and yr
    diagInd = logical(eye(numEdges,numEdges));
    totalsVec(:,:) = durTtlVec (grp,yr,1,:);
    totalsMat(:,:) = durTtlMat (grp,yr,:,:);
    totalsMat(diagInd) = 0;
    totalsMat(diagInd) = -sum(totalsMat,2);
    nonEmpty = find(totalsVec>0);
    transMat(nonEmpty,:) = bsxfun(@ldivide,...
        totalsVec(nonEmpty),totalsMat(nonEmpty,:));
    
    % set subintervals
    subinterval = bsxfun(@rdivide, transMat, -transMat(diagInd));
    [vs, vi] = sort(subinterval,2); % vs is sorted vectors, vi is indexed vectors
          
    % simulation
    display(['start simulation for group' num2str(grp)]);
    
    for i = 1: numSimulate/100
        disp( strcat(['Group' num2str(grp) ': i= ' num2str(i) ' x 100 interats']));
        for j = 1: 100
           
            for f = 1: length(exposureInd)
                
                t = 0;
                startR = ( cell2mat(exposureInd(f,2)) - 1 ) / 0.5;
                t = t+ exprnd(-1/transMat(startR,startR));
                
                while t < T
                    x = unifrnd(0,1);
                    jumpR = jumpRate(x, vs, vi, startR);
                    if ~(jumpR == 8)
                        startR = jumpR;
                        t = t+ exprnd(-1/transMat(startR,startR));
                    else
                        loss_portflio(100*(i-1)+j,grp) = loss_portflio(100*(i-1)+j,grp) + cell2mat(exposureInd(f,3))* (1-R);
                    end
                end   
            end
            
            
        end
    end
    
    loss = loss_portflio(find(loss_portflio)>0,grp);
    subplot(2,2,grp);
    hist(loss);
    
    
    disp(strcat(['fini group' num2str(grp)]));
    
end

figure
for grp = 1:4
loss = loss_temp(find(loss_temp)>0,grp);
[y,xi] = ksdensity(loss);
subplot(2,2,grp);plot(xi,y);
end


